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ABSTRACT 

This paper presents numerical techniques for constructing nonlinear predictive models to 
forecast solar flux directly from its time series. This approach makes it possible to extract 
dynamical invariants of our system without reference to any underlying solar physics. We 
consider the dynamical evolution of solar activity in a reconstructed phase space that 
captures the attractor (strange), give a procedure for constructing a predictor of future solar 
activity, and discuss extraction of dynamical invariants such as Lyapunov exponents and 
attractor dimension. 
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1. INTRODUCTION 
1 .1 Review of Solar Activity 

Need for Solar Flux Prediction. Solar flux F 10 .7 [radio flux emitted at a wavelength of X = J0.7cenume- 
ters (cm)] is the best indicator of the strength of ionizing radiations, such as solar ultraviolet and X- y 
emissions, that directly affect the atmospheric density and thereby change the orbit lifetime of satelhtes. 
Thus, accurate forecasting of solar flux F 10.7 is crucial for orbit determination of spacecrafts. 

Sunspots and Solar Flux. The strong correlation between sunspots and the solar flux F 10 . 7 is probably due to 
the enhanced radiation from limited areas of the Sun where sunspots are active. Sunspot activity depcndson 
the wavelength of radiated solar flux. For waves shorter than 3 cm, the intensity is steady. From 3 to 60 cm, 
often calledtoe decimeter range, the intensity occasionally increases for a few mmutes^Rising from the 
vicinity of active sunspot regions, decimeter intensity also tends to exhibit a 27-day period associated th 
solar rotation (Reference 1). 

The dynamics of sunspots and their formation are still a mystery. They are often more than 1000 degrees 
Kelvin cooler than the surrounding photosphere. Although many explanations for sunspot cooling have been 
proposed (the Biermann field inhibition mechanism and the superadiabatic downflow mechanism), the uge 
difference in temperature between sunspots and their surroundings suggests a similarity with solitons of 
multilevel turbulence. One may think of sunspots as solitons in a Tluid turbulent 2 ^° blt 

lifetime is a function of atmospheric drag, which is a function of atmospheric denstyw^chmt^isa 
function of solar flux. For this reason, spacecraft orbit determination requires accurate forecasting of solar 

flux. 


* This work was supported by the National Aeronautics and Space Administration (NASA)/Goddard Space Right Center 
(GSFC), Greenbelt, Maryland, Contract NAS 5-31500. 
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Nonlinear Structure in Solar Flux. Until recently, we had little reason to doubt that weather is in principle 
predictable, given enough data. Recently, a striking discovery changed our perspective: simple deterministic 
systems with only a few degrees of freedom can generate random behavior. When a system exhibits apparent 
random behavior that is fundamental to its dynamics, such that no amount of information gathering will make 
the system predictable, the system is considered to be chaotic. Much evidence supports our assertion that solar 
flux signal falls in this category (References 3 through 11). Perhaps paradoxically, chaos is generated by fixed 
rules that do not themselves involve any element of chance. Theoretically, the future of a dynamic system is 
completely determined by present and past conditions. In practice, however, amplification of small initial 
uncertainties makes a system with short-term predictability unpredictable in the long term. 

Many people speak of random processes as though they were a fundamental source of randomness. This idea 
is misleading. The theory of random processes is an empirical method to deal with incomplete information- it 
does not attempt to explain randomness. As far as we know, the only truly fundamental source of randomness 
is the uncertainty principle of quantum mechanics;. everything else is deterministic, at least in principle. 
Nonetheless, we call many phenomena, such as solar dynamics, random, even though we may not ordinarily 
think of them in terms of quantum mechanics. Historically, scientists have assumed that randomness derives 
solely from complication. In this paper, we will take the practical position that randomness occurs to the 
extent that a system’s behavior is unpredictable. We believe that randomness is subjective and a matter of 
degree; that is, some systems are more predictable than others (e.g., solar activity is more predictable than 
geomagnetic activity). 

Solar Activity Prediction. Interest in solar activity has grown in the past two decades for many reasons. 
Some reports claim a correlation between solar activity and weather on Earth (Reference 12), although a 
correlation has not yet been convincingly established (Reference 13). We have some evidence for the 
c oinc ident occurrences of the Maunder minimum (a period of little or no solar activity occurring from 1 645 to 
1715) and the “Little Ice Age” (a period of abnormally cold weather) (Reference 14). Perhaps most 
importantly for flight dynamics, solar activity changes the atmospheric density, which has important 
implications for spacecraft trajectory and lifetime prediction (Reference 15). The see min gly random nature 
of solar flux has misled us for many decades, causing us to assume that the underlying physics must 
necessarily be complex as well. Therefore, researchers have generally used statistical models to predict solar 
activity (Reference 16). However, new developments in chaos and nonlinear dy nami cs allow us to model the 
behavior of a system in terms of some invariants directly extractable from system dynamics, without 
reference to any underlying physics. Using chaos theory, we can predict short-term activity more accurately 
than with statistical methods; however, chaos theory imposes a fundamental limit on long-term predictions. 

1 .2 Brief Review of Chaotic Dynamics 

Self-Organization and Attractors. Imagine a very simple system: a pendulum. The pendulum exhibits two 
fu nd a m ental degrees of freedom: position and momentum. However, in its stable periodic state ( limi t cycle), 
the pendulum can be described by only one degree of freedom, the phase angle. Here, the dynamic is attracted 
to a lower-dimensional phase space, and the dimension of this reduced phase space is equal to the number of 
active degrees of freedom in the self-organized system. 

Attractors are not limited to zero dimension (fixed point) or one dimension (limit cycle), but for nonlinear 
systems they could be high dimensional and in some cases even fractional or fractal (strange attractors). 

Nonlinear Dynamical Systems. Anything that moves or evolves in time is a dynamical system. (If it does not 
move, it is a dynamical system at a fixed point.) Mathematically speaking, a dynamical system can be 
represented by a state space (phase space) R M and an evolution operator \[/ that defines how the state of the 
system evolves in time. M is the number of degrees of freedom in the dynamics; can be visualized as the 
physics of the system. 
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The attractor A of a dynamical system is the subset of phase space toward which the system evolves: 

'W(S) = A. 


lim 

t -*> 00 


An initial condition X 0 that is sufficiently near the attractor will evolve in time so that y, (XJ gets very close 
to the set A as t — ► oo . 

Phase-Space Construction Directly From a Time Series. When confronted with a complicated physical 
system, an experimenter normally measures at regular and discrete intervals of tone the value of some state 
variable (e.g., flux Fio 7 ) and records the time series f(to), f(ti), ffo),-. withf(ti) G R and ti to+iAt. Fro 
observed time series, the experimenter attempts to infer something about the dynamics (i.e., the physics) of 
the system. The measure f(t) represents a projection tc from R to R: 

jt: R m -+ R. 

Because the tone series is one-dimensional, it is an incomplete description of a system during a time 
evolution Nonetheless, many features of the dynamics can still be inferred from the tone senes alone, 
time-delayed values of the scalar time series, Takens (Reference 17) and Packard et al. (Reference 8) have 
shown that one can embed the time series into a higher dimensional space. Vectors are created wi 

components as 

f(t) = [f(t), f(t - t), ... f(t - (m - 1 )t)] t , 

wherex (time delay) and m (the embedding dimension) are parameters of the embedding procedure. Here f (t) 
represents a more complete description of dynamics than f(t) and can be thought of as a mapping. 


71 


m :R m — f R m . 


An embedding dimension of m > 2D+1, where D is the fractal dimension of the attractor, almost always 
ensures the construction of the topology of the attractor (Takens’ theorem. Reference 17). 

If unlimited infinitely precise data are available, almost any delay time x and embedding dimension m>D 
w iH w ork, atleast in principle. However, choosing the optimal parameters for real data is a nontnvialprocess. 

For example if x is too large then the components f(t) and f(R(m-l)x) of the reconstructed vector f will be 
effectively ^correlated, which will inflate the estimated dimension. On the other hand, if (m-l)x is too > small, 
then the components f(t), ...f(t+(m-l)x) will all be very nearly equal, and the reconstructed attractor will loo 
like one longdiagonal line. Generally, x must not be less than some characteristic decorrelation time, and 
(m-1 Yi musfnot be much greater than this decorrelation time. One such characteristic tone is the local minima 
of theautocorrelation function R(x) = <(f (tHOX^HO)). where < > represents average over tone. 

1 .3 Some Invariants of Dynamical Systems 

Lyapunov Exponent In a chaotic system, the adjacent points of the time series become separated under the 
action of a map; in our case, f n is the value of solar flux measured daily. 

fn+1 = M(fn)> 

Which leads to satic motion. The Lyapunov exponent measures this exponential separation, as shown in 

Figure 1. 

Therefore, 

e 0 e^W = |M N (f 0 + e 0 ) " mN ( f o)l- 
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Figure 1. Exponential Separation Measurement 
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0 and N -* oo , we get the formal expression for Lyapunov exponent Xff 0 ): 
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The Lyapunov exponent also measures the average loss of information. 

In variant Measure. The invariant measure Q(f) determines the density of the iterates of a map over the unit 
interval and is defined to be 


e(t) 


= i?-sl 5 [ f - M’ow], 

i = 0 


Kolmogorov Entropy. Kolmogorov entropy K describes the dynamical behavior at the strange attractor K is 
the analog of thermodynamic entropy that measures the disorder in a dynamical system. For a 
one-dimensional map, it is just the Lyapunov exponent. The rate K (at which information about the system is 
lost) is equal to the average sum of positive Lyapunov exponents: 


k = J (f), 

+ i 

where the superscript d is the dimensionality of our phase space. 

In most cases, the As are independent of f, so 


K -Zvw- 


K is, indeed, a useful measure of chaos. K becomes zero for regular motion. It is infinite for totally stochastic 
or random systems but is a constant larger than zero for chaotic systems (Figure 2). In higher dimensional 
systems, we lose information about the system (as shown in Figure 3). Here, a sphere of radius £o changes its 
geometry in phase to an ellipsoid as the system evolves in time. 


Ha usd orff Dimension. One of the invariants of an attractor that can be extracted directly from the time series 
is called Hausdorff Dimension D, an infinite set of dimensions D=D 0 , Dj, D 2 , ... that describes the 
omogeneity of the attractor. It can be shown that D 2 (which yields a lower bound on the Hausdorff 
dimension) and many other invariants of the system can be directly obtained from a time series: 

• D 2 , that is, a lower bound on the Hausdorff dimension (D 2 < D]) 

• d, that is, the embedding dimension of the attractor 

• The amplitude of the white noise on the signal; that is, separating the deterministic chaotic 
motion of the system from disturbing white noise 


322 




Figure 2. Kolmogorov Entropies tor Regular, Chaotic, and 
Random Motions 
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Figure 3. A Two-Dimensional Map of Small Circle Into an Ellipse 
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A lower bound on the Kolmogorov entropy; that is, determining “how chaotic” the signal is 

* tension: connects the static structure of the attractor (as measured by 

U-D 0 , D, , D 2 , ...) and the dynamical behavior at the attractor (as measured by the Lyapunov 
exponents) 


Dimensions of a Strange Attractor. To characterize the inhomogeneous static structure of the 
introduce an infinite set of dimensions D n related to the nth powers of Pj via 


attractor, we 


log 


— £, 


lim 


Q(b 0 ) 

1 i-0 


n - 1 


loge c 


n = 0, 1, 2, 

where Pi is the probability of finding a point of attractor in the cell number i [i=l , 2, ... Q (e 0 )l . For n=0 we vet 
e Hausdorff dimension of the attractor D=D 0 ; for n=l, we get dimension D=Di (called the information 

dimension). It should be noted that for n -+oo , the dimension is still a finite number (Definite) which in 
general is not an mteger. t% 

Largest Lyapunov Exponent of Solar Flux Time Series. The sum of the Lyapunov exponents is the 
ime-averaged divergence of the phase space trajectory; hence, any dissipative dynamical system will have at 
least one negative exponent. Any dynamical system without a fixed point will have at least one zero Lyapunov 

A small positive :Lyapunov exponent is an indication ofehaos, and a very large positive Lyapunov exponent is 

an mdicauon of a totally stochastic or random system. Therefore, the sign of the exponent prides a 
qualitative picture of a system’s dynamics— a positive exponent represents chaos, a zero exponent represents 

S?? ? / ySte “ S ’ “ d 3 negative ex P° nent ^presents periodic systems. Figures 4 and 5 show the 
actual solar flux data and the largest Lyapunov exponent, respectively, for more than 4000 points. Here we 
have used the well known technique of phase space reconstruction with delay coordinates (Reference 18) 
Mter embedding the solar flux tune series in a state space using the Takens-Packard delay coordinate 
echmque, one can learn the induced nonlinear mappmg using a local approximation. This will allow us to 
make short-term forecasting of the future behavior of our time series using infoimation based only on past 

(Refe S The i55 C ° r eSdmate 0f such a techni <l ue h as already been developed by Farmer and Sidorowich 


E = Ce( m+1 ) KT N-(m+t)/D 

9 


where E 

m 

K 

T 

N 

D 

C 


normalized error of prediction (0$E< 1, where zero is perfect prediction and one 
is a prediction no better than average) 

order of local approximation 

Kolmogorov entropy 

forecasting window 

number of data points 

dimension of the attractor 

normalization constant 


Usmg the Farmer-Sidorowich relation, we can find the prediction horizon T for the zeroth order of local 
approximation. Any prediction above T max is no better than average constant prediction. 
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Figure 4. Actual Solar Flux Time Series 
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Figure 5. Largest Lyapunov Exponent of Solar Flux Time Series of Figure 4 


E(T m ax) = 1 

Thus, for m = 0, K is the largest Lyapunov exponent X. Therefore, 

, m In (N) 

gKTmu J\J l/D _ 1 or T max 

MN) 

lmax XD 


325 




Any prediction beyond the indicated horizons is no better than average value. The connection between the 
local and the global Lyapunov exponents has recently been found (March 1991) by Abarbanel et al. 
(Reference 20) in a form of power law as 

M 1 ) = + p 


N = col. 


where X(l) 
1 

v 

c 

CO 


— local Lyapunov exponent 
“ length of observed data (observation window) 

= a constant dependent to the dynamical system (0.5 <v< 1.0) 
= a constant dependent to initial conditions of the system 
= well known global Lyapunov exponent 
= frequency of data points 


Because any data are finite length data, using the Abarbanel-Kennel power law and Farmer-Sidorowich 
relation, we can find T max as 


T 


max 


ln(lco) 

(*° + f) d ' 


This means that as 1 increases linearly, T max increases logarithmically to a certain asymptotic T because of the 
denominator c/l v (Figure 6). 

Therefore, our relation shows that at the asymptote T raax = T 0 (Reference 1 1 ) and dT max /dl = 0. Thus, we can 
fmd what observation window is required for forecasting up to T max within some confidence level. 
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dN 


* 0 («> 

0, thus N 0 ~ e~ , x 0 (5) > 2, 



* max \ 

/ 

'o 

1 Q max \ 

/ : 

/ 

t, 

i 

L ! 

/ 

. '« i 

/ “ — 

/ 

T 3 ! 

r 

/ 

' 3 

/ 

'o 


= finite 


Figure 6. Relationship of I to T max 
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where xq (5) is the solution to e x (x - 1) = 5, and where 



is the scaled global Lyapunov exponent. 

This result shows that any observation window greater than lo = No/co will not improve our predic tion horizon 
T 0 ; so more data beyond this limit are not needed to understand a dynamical system. This conclusion is mdeed 
consistent with weather prediction and also with empirical results concluded from neural networks training. 

2. STRUCTURE OF THE COMPUTER PROGRAM 

Once we know the state space representation, the next goal is to fit a model to the data. There are several 
methods. The simplest method is to assume that the dynamics can be written as a map in the form 

f n+1 - M(f n ), 

where the current state is f n , and f n+ i is a future state. Methods such as the polynomial method, rational 
approximation, radial basis functions, neural networks, and local approximations have been proved to be 
successful. Here we only introduce local approximation technique, which is the method used to structure the 

computer program. 

Local Approximation. The basic idea is to break up the domain of M into local neighborhoods and fit 
different parameters into each neighborhood. This fit is generally better than global approximation 
especially for large data sets. Most global representations reach a point of diminishing returns, at which 
adding more parameters or data gives only an incremental improvement in accuracy. After a certain point, 
adding more local neighborhoods is usually more efficient than adding more parameters and going to higher 
order. With local approximation, it is possible to use a given functional representation efficiently. The key is 
to choose the local neighborhood size correctly, so that each neighborhood has just enough points to make the 
local parameter fits stable. The basic idea is shown in Figure 7 . 

Moving to representations of higher degree involves a tradeoff— higher degree representations promise more 
accuracy, but also require larger neighborhoods. A larger neighborhood implies that the complexity of M 
increases. Finding the best compromise between these two effects is a central issue in local approximahon. 

An example of local approximation is first order, or nearest neighbor, approximation. Look through the data 
set for nearest neighbor, and predict the current state based on what the neighbor did at time T later. We 
approximate f(t + T) by f(t, T) = f(t' + T), where f(f), is the nearest neighbor of f(t). That is, to prechct 
tomorrow’s solar flux, we would search the historical record and find the solar flux pattern most similar to that 
of today, and predict that tomorrow’s solar flux pattern will be the same as the neighboring pattern 1 day later. 
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First order approximation can sometimes be improved by finding more neighbors and merging their 
predictions, for example, by weighting according to distance from the current state. If the data are noise, it is 
better to use a larger number of neighbors. This procedure can be improved by weighting the contributions of 
neighbors according to their distance from the current state. The beauty of linear approximation is that the 
neighborhood size grows slowly with the embedding dimension. The order of approximation may depend on 
factors such as the choice of neighborhoods, the dimension, or peculiarities of the data set. For low 
dimensional problems, a third order approximation is good. 

Nonstationary Behavior. If the trajectory is on an attractor, the data must be stationary, as long as the 
parameters are constant. However, parameter variations can result in nonstationary behavior. To deal with this 
problem, we include time as one of the state space coordinates. Time can be included in the metric, so that the 
search for nearest neighbors favors recent data. This approach takes into account the trends, other 
time-dependent effects, and seasonality. 

Implementation of Local Approximation. Finding neighbors in a multidimensional data set is time 
consuming when considering many points. A straightforward method is to compute the distance to each 
point, which takes approximately N steps for N points. This can be reduced to roughly log N steps by 
organizing the data with a decision tree, such as a k-d tree (Reference 21). 

In this method, die data set is partitioned one coordinate at a time. We can take the coordinate with largest 
range and partition it at its median value. These values are stored in the tree as keys. It is now possible to 
eliminate many points form consideration when looking for the nearest neighbors. This way, we minim i 
processing time considerably. 

2.1 CHOICE OF THE EMBEDDING DIMENSION d 

Here we would like to determine the correct value of the embedding dimension d from the scalar time series 
x(n), n=l ,2,...,Nd. We assume that there are enough data that we need not be concerned with statistical issues 
about numerical accuracy. We also assume that extrinsic noise is absent from the x(n) when we receive them. 
We further assume that by following Takens’ phase-space reconstruction technique we have successfully 
captured the dynamics and embedded our time series. This requires a correct choice of t, which will be 
discussed in the next section. 

For now, let’s further assume we have a correct t to construct the attractor in the phase space. To establish 
dimension d, we need some characteristic of the attractor that becomes unchanging as d becomes large 
enough, thus indicating that the attractor can be embedded in R d . This invariant characteristic of the attractor 
is the attractor dimension d A . One increases d until d A remains constant and identifies the minimum d where 
d A “saturates” as the embedding dimension. But computation of d A is difficult, so we use the correlation 
function D(r) proposed by Takens (Reference 17). 

D(r,N ' d) = N(N 2 - 1) I I U < r HI*© " f(0 II) i * j, 

where U(z) is just the unit step function U(z) = J . 

For N large enough, the behavior of D(r,N,d) for r becomes independent of N and D(r,N,d) takes the form 

D(r, N, d) = d>(r,d)r v < d ). 

If we plot D(r,N,d) versus r we can single out the correct value of dimension d as in Figure 8. From Figure 8 it 
is concluded that the minimum value of d=3 is the right choice beyond which attractor dimension d A does not 
change or the slope of our graph becomes constant. 

In the next section we study the correct choice of x to reconstruct the phase-space attractor. 
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Figure 8. Correlation and Attractor Dimensions 
2.2 Choice of the Time Shift X 

The choice of time shifts x is not well agreed upon. If the underlying system is described by a differential 
equation and scalar variable x(t) were measured at discrete times x(n) =x(to + n At), then we would be, by the 
choice of lagged variables, trying to find a discrete replacement for the usual phase-space coordinates: 

x(t),x(t), 

The best choice for time shift t is a fraction of the time associated with the first local minimum of the 
autocorrelation function 

^ | x(t + x)x(t)dt . 

We find that this choice, although somewhat arbitrary, works well in practice and gives a simple systematic 
way of determining x. 

3. RESULTS 

Figure 9 shows the daily observed values of solar flux F 10.7 for about four solar cycles from Februaiy 1 947 to 
November 1991. A close examination of this graph shows low daily variability at solar cycle minima and 
large daily variability at solar cycle maxima. Therefore, the challenge for solar physicists is to forecast solar 
flux in the region of solar cycle maxima. 

Here forecasts are made in the region of solar cycle maxima, where the variations are as large as 1 50 units of 
solar flux and in regions between a maximum and a minimum, where the variations are as large as 100 units 
of solar flux. Comparisons are made versus 27-day NOAA predictions, Schatten’s monthly predictions, and 
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observed values of flux F 10 . 7 . Figure 1 0 shows 27 -day predictions using chaos theory with undecoupled (raw) 
data. Predictions were made in the region of solar cycle maxima. Comparisons show clearly that chaos 
predictions are at worst 20 units off from observed values, whereas the best NOAA 27-day predictions are 
about 40 units off from observed values, with a wrong phase on the 27-day local maximum. 

Here it is clear that chaotic long-term prediction is very close to the 81-day average of actual Fjoj, which is 
useful for calculating exospheric temperature using the Jacchia-Roberts atmospheric density model. 

Figure 1 1 shows the same conclusion for predictions made for a 60-day span. Here it is very clear that NOAA 
predictions (best updated values to the end of 60 days) do not cany any of the significant dynamical features 
of the observed values, whereas the chaos prediction does carry the inherent dynamical features. Figure 12 
shows a comparison of MSFC and Schatten’s predictions. Figure 13 also shows comparisons of 81 -day 
average Fjo. 7 , Schatten’s predictions, and chaos predictions. 

Figure 14, which was made for a region between solar cycle maxima and minima, shows a 30-day prediction 
indicating that chaos prediction is at worst 20 units off from the observed values, but it should be noted that 
the signal-to-noise ratio characteristic of solar activity indicators are - 35 forF] 0 . 7 , - 20 for sunspot numbers, 
and - 10 for sunspot areas. 



Figure 9. Daily Observed Values of Solar Flux F 10 .7 (Watt per m 2 per Hz x 
10 -22 ) Versus Modified Julian Date 



Figure 10. 27-Day Chaos and NOAA Predictions Compared With Observed 
Solar Flux F 10 7 From April 13, 1990 
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Figure 11. 60-Day Chaos and NOAA Predictions (Best Weekly Updates) 
Compared With Observed Solar Flux F 10.7 
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Figure 12. Comparison of MSFC and Schatten Solar Flux F 10 .7 Predictions With 
Observed F 10 .7 
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Figure 13. 30-Day Comparison of 81 -Day Average of Observed F 10 .7> Schatten’s 
Prediction and Chaos Prediction With Observed Daily F 10 .7 








Figure 14. 30-Day Comparison of Schatten’s Prediction and Chaos Prediction 
With Observed Daily F 10 .7 

Figure 15 shows chaos predictions for embedding dimensions 3 and 4, the Schatten mean, and +2o 
predictions. Figure 15 clearly shows that, for long-term forecasting, the chaotic method carries all the 
inherent structures of Schatten’s method. Figure 16 shows that chaos prediction captures the cyclical behavior 
of solar cycles (the minimum of the cycle is very clear). Figure 17 shows the observed F 10 .7 and its 81 -day 
average, and Figure 18 shows chaos prediction for about 50 months after November 1991 (the time of our 
analyis). Comparsions of chaotic predictions (Figure 1 6) with Schatten’s predictions show the chaos model to 
be predicting F 10 .7 even higher than Schatten’s +2a. Recently, the observed F 10 .7 have in fact been higher than 
Schatten’s +2 o. As seen in Figure 17, after October the average flux is about 200, as is clear from the first 
couple of points in Figure 18. It is certainly possible to fine tune the model by adjusting the embedding 
dimension and the time shift t. Figure 19 shows predictions for various time shifts x, with D=3. Figure 20 
shows predictions for various time shifts x, with embedding dimension D=4. 

4. CONCLUSIONS 

In this paper we presented numerical techniques for constructing nonlinear predictive models to forecast solar 
flux F10.7 directly from its time series. Using this approach, we extracted dynamical invariants of our system 



Figure 15. Chaos Long-Term Prediction With Embedding Dimension D=3 and 
D=4 Compared With Schatten’s Mean and Schatten’s +2 a 


332 






Months From Nov 1 991 


Figure 16. Predictions of Schatten and Chaos Theory for Solar Cycle Minima 



Figure 17. Observed Solar Flux F 10 .7 and Its 81 -Day Average, Starting From 
October 1 991 



igure 1 8. 50-Month Chaos Prediction With Embedding Dimension D=4 and 
Time Shift x =10 Days 








Figure 19. Variability of Predictions as a Function of Change in Time Shift x. 

Here Embedding Dimension Is Kept Constant (D=3), and x Is Varied 
by 8, 10, 17, and 21 Days 



Figure 20. Variability of Predictions as a Function of Change in Time Shift x. 

Here Embedding Dimension Is Kept Constant (D=4), and x Is Varied 
by 8, 10, 17, and 21 Days 

without reference to any underlying solar physics, thereby circumventing the complicated physics and 
modeling the system directly from data. 

Comparison of our prediction of solar flux activity using chaotic dynamics with conventional methods used 
by Schatten, MSFC, and NOAA demonstrated the validity of our approach to modeling solar activity using 
nonlinear dynamics. This approach could also be used to model other complicated systems, such as 
geomagnetic activity and atmospheric density, to name just two. 
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